home *** CD-ROM | disk | FTP | other *** search
/ CU Amiga Super CD-ROM 19 / CU Amiga Magazine's Super CD-ROM 19 (1998)(EMAP Images)(GB)[!][issue 1998-02].iso / CUCD / Programming / LEDA / source / src / plane / _delaunay_tree.c < prev    next >
Encoding:
C/C++ Source or Header  |  1994-11-16  |  41.5 KB  |  1,558 lines
  1. /*******************************************************************************
  2. +
  3. +  LEDA  3.1c
  4. +
  5. +
  6. +  _delaunay_tree.c
  7. +
  8. +
  9. +  Copyright (c) 1994  by  Max-Planck-Institut fuer Informatik
  10. +  Im Stadtwald, 6600 Saarbruecken, FRG     
  11. +  All rights reserved.
  12. + 
  13. *******************************************************************************/
  14.  
  15.  
  16. /*****************************************************************************
  17. +
  18. +    Delaunay tree    by Olivier Devillers
  19. +
  20. +        Fully dynamic construction of the Delaunay triangulation with
  21. +        good randomized complexity
  22. +
  23. +    Copyright (c) 1990
  24. +    INRIA, 2004 route des Lucioles, 06 565 Valbonne Cedex, France
  25. +
  26. +
  27. ******************************************************************************/
  28.  
  29.  
  30. #include <LEDA/impl/delaunay_tree.h>
  31.  
  32.  
  33. #define Fini        0
  34. #define Infini1     1        /* 1 point a l'infini */
  35. #define Infini2        2        /* 2 points a l'infini */
  36. #define Infini3        3        /* 3 points a l'infini, c'est la racine */
  37. #define Impossible    99        /* Le voisin du precedent */
  38. #define U 0                /* indices pour les trois sommets */
  39. #define V 1
  40. #define W 2
  41.  
  42.  
  43. /*****************************************************************************
  44.  
  45.             Declarations de type
  46.  
  47. ******************************************************************************/
  48.  
  49.  
  50. struct LISTENOEUD{
  51.     NOEUD *item;
  52.     LISTENOEUD *next;
  53. };
  54.  
  55.  
  56.  
  57. struct REINSERER {
  58.     DT_item    x;                /* a point devant etre reinserer */
  59.     NOEUD   *detruit1,            /* triangles a detruire si ils existent */
  60.         *detruit2;            /* NULL sinon xpx1 est direct et xpx2 indirect*/
  61.     listenoeud
  62.         decroches;            /* triangles decroches */
  63.     REINSERER
  64.         *finf,*fsup;        /* pour l'arbre binaire */
  65. };
  66.  
  67.  
  68.  
  69. struct STAR{
  70.     STAR
  71.             *next,*previous;/*arete suivante et precedente sur le polygone*/
  72.     NOEUD
  73.             *triangle;        /* triangle adjacent */
  74.     int        n,p,arete;        /*point suivant, precedent, arete dans le triangle*/
  75. };
  76.  
  77.  
  78.  
  79.  
  80. struct NOEUD{
  81.  
  82.     int        type;    /* le triangle est Fini ou Infini[123] 
  83.                    dans le cas d'un triangle infini, les sommets
  84.                    infinis sont la direction a` l'infini */
  85.     
  86.     int        degenere;  /* booleen, 3 pts alignes ? */
  87.     int        direct;       /* booleen, le triangle est-il direct ? */
  88.         int             stacked;
  89.     int        visite;       /* triangle marque visite si ce flag
  90.                           est le flag courant.
  91.                       pour la suppression le triangle est
  92.                                   marque si il appartient ou a
  93.                                   appartenu a A 
  94.                                 */
  95.     
  96.     int        detruire;   /* le triangle est a detruire si ce flag
  97.                        est le flag courant */
  98.     coord     cx,cy;            /* centre du cercle 
  99.                                           (uniqu^t pour tracer voronoi)*/
  100.     coord    a2,b2,a3,b3,c1,z0,z1,z2; /* pour les tests de conflit */
  101.     
  102.     DT_item    meurtrier; /* NULL si le triangle n'est pas mort */
  103.     DT_item    s[3];       /* sommets du triangle */
  104.                    /* U est le createur */
  105.                    /* si Infini1 alors W est a l'infini */
  106.                    /* si Infini2 alors V et W sont a l'infini */
  107.                    /* si Infini3 alors U,V et W sont a l'infini*/
  108.     NOEUD
  109.             *pere,        /* pere */
  110.             *bpere,        /* beau-pere */
  111.             *fils[3];    /* fils */
  112.     listenoeud
  113.             bfils;        /* liste des beaux-fils */
  114.     NOEUD
  115.             *voisin[3],    /* 3 voisins (courant ou a la mort) */
  116.             *creation[3],    /* 3 voisins a la creation */
  117.             *special[3];    /* 3 voisins speciaux! */
  118.     star    Star[3];        /* pointeur reciproque vers Star */
  119.  
  120.  
  121. };
  122.  
  123. /*****************************************************************************
  124.  
  125.             Allocation memoire
  126.  
  127. ******************************************************************************/
  128.  
  129. #define SIZEALLOC 100
  130.  
  131.  
  132. char* delaunay_tree::alloc_block(int size)
  133. {
  134.    char* p =  (char*)malloc( size + sizeof(double) ); 
  135.  
  136.    if (!p) error_handler(1,"out of memory");
  137.  
  138.    *((char**)p) = used_blocks;
  139.    used_blocks = p;
  140.  
  141.    return (used_blocks+sizeof(double)) ; 
  142.  }
  143.  
  144.  
  145. void delaunay_tree::free_blocks()
  146. {
  147.   char* p;
  148.   while (used_blocks)
  149.   { p = used_blocks;
  150.     used_blocks = *((char**)used_blocks);
  151.     free(p);
  152.    }
  153. }
  154.  
  155. void delaunay_tree::poubelle_point(DT_item p)
  156. { p->cree = (noeud) Poubelle_point;
  157.   Poubelle_point = p;
  158. }
  159.  
  160. DT_item delaunay_tree::vrai_alloc_point()
  161. {
  162.  if ( ! item_nombre ) {
  163.         item_next = (DT_item) alloc_block( (item_nombre = SIZEALLOC) * sizeof(DT_POINT) ) ; 
  164.     }
  165.     item_nombre--;
  166.     return item_next++;
  167. }
  168.  
  169. DT_item delaunay_tree::alloc_point()
  170. {
  171.     if ( Poubelle_point ) { DT_item p = Poubelle_point;
  172.                 Poubelle_point = (DT_item) Poubelle_point->cree;
  173.                 return p;  }
  174.     return vrai_alloc_point();
  175. }
  176.  
  177.  
  178. void delaunay_tree::poubelle_listenoeud(listenoeud p)
  179. {
  180.   p->next =  Poubelle_listenoeud;
  181.   Poubelle_listenoeud = p;
  182. }
  183.  
  184. listenoeud delaunay_tree::vrai_alloc_listenoeud()
  185. {
  186.  if ( ! listenoeud_nombre ) {
  187.         listenoeud_next = (listenoeud) alloc_block((listenoeud_nombre = SIZEALLOC) * sizeof(LISTENOEUD) ) ; 
  188.     }
  189.     listenoeud_nombre--;
  190.     return listenoeud_next++;
  191. }
  192.  
  193. listenoeud delaunay_tree::alloc_listenoeud()
  194. {
  195.     if ( Poubelle_listenoeud ) 
  196.                               { listenoeud p = Poubelle_listenoeud;
  197.                     Poubelle_listenoeud = Poubelle_listenoeud->next;
  198.                 return p;  }
  199.     return vrai_alloc_listenoeud();
  200. }
  201.  
  202.  
  203. void delaunay_tree::poubelle_noeud(noeud p)
  204. /* on ne se resert que des noeud dont le champs detruire n'est pas flag */
  205. { 
  206.   p->fils[0] =  Poubelle_noeud;
  207.   Poubelle_noeud = p;
  208.   if (fl != p->detruire ) { fl = p->detruire; Libre_noeud = Poubelle_noeud ;}
  209. }
  210.  
  211. noeud delaunay_tree::vrai_alloc_noeud()
  212. {
  213.  if ( ! noeud_nombre ) {
  214.         noeud_next = (noeud) alloc_block((noeud_nombre = SIZEALLOC) * sizeof(NOEUD) ) ; 
  215.     }
  216.     noeud_nombre--;
  217.     return noeud_next++;
  218. }
  219.  
  220. noeud delaunay_tree::alloc_noeud()
  221. /* on ne se resert que des noeud dont le champs detruire n'est pas flag */
  222. {    
  223.     if (Libre_noeud && Libre_noeud->fils[0] )
  224.             { noeud p = Libre_noeud->fils[0];
  225.               Libre_noeud->fils[0] = Libre_noeud->fils[0]->fils[0];
  226.               return p;  }
  227.     return vrai_alloc_noeud();
  228. }
  229.  
  230.  
  231. void delaunay_tree::poubelle_star(star p)
  232. {
  233.   p->next =  Poubelle_star;
  234.   Poubelle_star = p;
  235. }
  236.  
  237. star delaunay_tree::vrai_alloc_star()
  238. {
  239.  if ( ! star_nombre ) {
  240.         star_next = (star) alloc_block((unsigned int) (star_nombre = SIZEALLOC) * sizeof(STAR) ) ; 
  241.     }
  242.     star_nombre--;
  243.     return star_next++;
  244. }
  245.  
  246.  
  247. star delaunay_tree::alloc_star()
  248. {
  249.     if ( Poubelle_star ) { star p = Poubelle_star;
  250.                    Poubelle_star = Poubelle_star->next;
  251.                    return p;  }
  252.     return vrai_alloc_star();
  253. }
  254.  
  255.  
  256. void delaunay_tree::poubelle_reinserer( reinserer p)
  257. {
  258.   p->finf =  Poubelle_reinserer;
  259.   Poubelle_reinserer = p;
  260. }
  261.  
  262. reinserer delaunay_tree::vrai_alloc_reinserer()
  263. {
  264.  if ( ! reinserer_nombre ) {
  265.         reinserer_next = (reinserer) alloc_block((reinserer_nombre = SIZEALLOC) * sizeof(REINSERER) ) ; 
  266.     }
  267.     reinserer_nombre--;
  268.     return reinserer_next++;
  269. }
  270.  
  271.  
  272. reinserer delaunay_tree::alloc_reinserer()
  273. {
  274.     if ( Poubelle_reinserer ) {reinserer p = Poubelle_reinserer;
  275.                       Poubelle_reinserer = Poubelle_reinserer->finf;
  276.                    return p;  }
  277.  
  278.     return vrai_alloc_reinserer();
  279. }
  280.  
  281. /*****************************************************************************
  282.  
  283.             Structure Reinsert
  284.  
  285. ******************************************************************************/
  286.  
  287. reinserer delaunay_tree::locate( reinserer n, DT_item x)
  288. {
  289.     reinserer nn;
  290.  
  291.     if ( ! Reinsert ) {
  292.         nn = alloc_reinserer();
  293.         nn->x = x;
  294.         nn->decroches = NULL;
  295.         nn->detruit1 = NULL;
  296.         nn->detruit2 = NULL;
  297.         nn->finf = NULL;
  298.         nn->fsup = NULL;
  299.         Reinsert = nn;
  300.      return nn;
  301.     }
  302.     if ( n->x == x) return n;
  303.     if ( x->n < n->x->n ) {
  304.         if (n->finf) return locate(n->finf,x);
  305.         nn = alloc_reinserer();
  306.         nn->x = x;
  307.         nn->decroches = NULL;
  308.         nn->detruit1 = NULL;
  309.         nn->detruit2 = NULL;
  310.         nn->finf = NULL;
  311.         nn->fsup = NULL;
  312.         n->finf = nn;
  313.     } else {
  314.         if (n->fsup) return locate(n->fsup,x);
  315.         nn = alloc_reinserer();
  316.         nn->x = x;
  317.         nn->decroches = NULL;
  318.         nn->detruit1 = NULL;
  319.         nn->detruit2 = NULL;
  320.         nn->finf = NULL;
  321.         nn->fsup = NULL;
  322.         n->fsup = nn;
  323.     }
  324.     return nn;
  325. }
  326.  
  327. void delaunay_tree::clear_Reinsert( reinserer n)
  328. {    listenoeud l,ll;
  329.     if ( ! Reinsert ) return;
  330.     if (n->finf) clear_Reinsert(n->finf);
  331.     if (n->fsup) clear_Reinsert(n->fsup);
  332.     for ( l = n->decroches; l; l = ll ) {ll = l->next; poubelle_listenoeud(l);}
  333.     poubelle_reinserer(n);
  334.     if ( n == Reinsert ) Reinsert = NULL;
  335. }
  336.  
  337. /*****************************************************************************
  338.  
  339.             Structure Star
  340.  
  341. ******************************************************************************/
  342.  
  343. void delaunay_tree::clear_Star()
  344. {
  345.     star s;
  346.  
  347.     if ( ! Star ) return;
  348.     for ( s = Star->next ; s != Star ;)
  349.                 {s = s->next; poubelle_star(s->previous);}
  350.     poubelle_star(Star);
  351.     Star = NULL;
  352. }
  353.  
  354. void delaunay_tree::init_Star( noeud n, int i, int j, int k)
  355. /* l'arete i,j de n est sur le bord de star */
  356. {
  357.     star s;
  358.  
  359.     if (! n) {
  360.         Star->next = first_star;
  361.         first_star->previous = Star;
  362.         first_star = NULL;
  363.         return;
  364.     }
  365.     s = alloc_star();
  366.     s->triangle = n;
  367.     s->p = i;
  368.     s->n = j;
  369.     s->arete = k;
  370.     n->Star[k] = s ;
  371.     s->previous = Star;
  372.     if (! Star ) first_star = s;
  373.     else         Star->next = s;
  374.     Star = s;
  375. }
  376.  
  377. star delaunay_tree::loc_Star( DT_item xx)
  378. /* trouve l'arete telle que xx soit le point precedent */
  379. {
  380.     star s;
  381.  
  382.     s = Star;
  383.     while ( s->triangle->s[s->p] != xx ) s = s->next ;
  384.     return s;
  385. }
  386.  
  387. void delaunay_tree::split_Star( star n1,star n2)
  388. /* raccroche les n1 a n2 en viarant la partie intermediare */
  389. {
  390.     star s,ss;
  391.  
  392.     if (n1==n2) {
  393.         s = alloc_star();
  394.         s->previous = n1;
  395.         s->next = n1->next;
  396.         n1->next = s;
  397.         s->next->previous = s;
  398.         return;
  399.     }
  400.     for ( s = n1->next; s != n2 ; ){ss=s->next; poubelle_star(s); s=ss; }
  401.     n1->next = n2;
  402.     Star = n2->previous = n1;
  403. }
  404.  
  405. /*****************************************************************************
  406.  
  407.             Procedures
  408.  
  409. ******************************************************************************/
  410.  
  411. void delaunay_tree::beaufils( noeud bf, noeud bp)
  412. /* ajoute bf a la liste des beaufils de bp */
  413. {
  414.     listenoeud nov;
  415.  
  416.     nov = alloc_listenoeud();
  417.     nov->next = bp->bfils;
  418.     nov->item = bf;
  419.     bp->bfils = nov;
  420. }
  421.  
  422. void delaunay_tree::plusbeaufils( noeud bf, noeud bp)
  423. /* enleve bf a la liste des beaufils de bp */
  424. {
  425.     listenoeud n,nn;
  426.  
  427.     if ( bp->bfils->item == bf ) {
  428.         nn = bp->bfils;
  429.         bp->bfils = bp->bfils->next ;
  430.         poubelle_listenoeud(nn);
  431.         return;
  432.     }
  433.     for ( n = bp->bfils; n->next->item != bf; n = n->next );
  434.     nn = n->next ;
  435.     n->next = n->next->next ;
  436.     poubelle_listenoeud(nn);
  437. }
  438.  
  439. void delaunay_tree::demiplan( noeud t)
  440. /* calcul de la normale */
  441. {
  442.     t->type = Infini1;
  443.     t->degenere = 0;
  444.     if ( t->direct ) {
  445.             t->cx= t->s[U]->y - t->s[V]->y ;
  446.             t->cy= t->s[V]->x - t->s[U]->x ;
  447.     } else {
  448.             t->cx= t->s[V]->y - t->s[U]->y ;
  449.             t->cy= t->s[U]->x - t->s[V]->x ;
  450.     }
  451. }
  452.  
  453. void delaunay_tree::cercle( noeud t)
  454. {
  455.     register double a,b,c,d;
  456. /* calcul pour les tests de conflits */
  457.        t->a2 = t->s[1]->x - t->s[0]->x ;
  458.        t->b2 = t->s[1]->y - t->s[0]->y ;
  459.        t->a3 = t->s[2]->x - t->s[0]->x ;
  460.        t->b3 = t->s[2]->y - t->s[0]->y ;
  461.  
  462.        t->c1 = t->a2 * t->b3 - t->a3 * t->b2 ;
  463.        t->z1 = t->a2 * (t->s[1]->x + t->s[0]->x) + t->b2 * (t->s[1]->y+t->s[0]->y);
  464.        t->z2 = t->a3 * (t->s[2]->x + t->s[0]->x) + t->b3 * (t->s[2]->y+t->s[0]->y);
  465. /* calcul du centre et du rayon du cercle  */
  466.     a = t->s[U]->x * t->s[U]->x + t->s[U]->y * t->s[U]->y ;
  467.     b = t->s[V]->x * t->s[V]->x + t->s[V]->y * t->s[V]->y - a ;
  468.     c = t->s[W]->x * t->s[W]->x + t->s[W]->y * t->s[W]->y - a ;
  469.     d = 2*( (t->s[V]->x - t->s[U]->x) * (t->s[W]->y - t->s[U]->y)
  470.                 - (t->s[W]->x - t->s[U]->x) * (t->s[V]->y - t->s[U]->y));
  471.     if ( (t->degenere = (d==0) ) ) { 
  472.                         t->cx= t->pere->cx ;
  473.                         t->cy= t->pere->cy ;
  474.                     return;
  475.                 }
  476.     t->cx=(b*(t->s[W]->y - t->s[U]->y)-c*(t->s[V]->y - t->s[U]->y))/d;
  477.     t->cy=(c*(t->s[V]->x - t->s[U]->x)-b*(t->s[W]->x - t->s[U]->x))/d;
  478. }
  479.  
  480. int delaunay_tree::confondu( DT_item a, DT_item b)
  481. {
  482.     return ((a->x == b->x)&&(a->y==b->y));
  483. }
  484.  
  485. int delaunay_tree::direct( noeud n)
  486. /* teste si un triangle est direct */
  487. {
  488.     switch(n->type){
  489.     case Fini       : 
  490.         return ((n->s[V]->x - n->s[U]->x)*(n->s[W]->y - n->s[U]->y)
  491.                 - (n->s[V]->y - n->s[U]->y)*(n->s[W]->x - n->s[U]->x) > 0);
  492.     case Infini1     :
  493.         return ((n->s[V]->x - n->s[U]->x)*n->s[W]->y
  494.                 - (n->s[V]->y - n->s[U]->y)*n->s[W]->x > 0);
  495.     case Infini2     :
  496.         return (n->s[V]->x*n->s[W]->y - n->s[V]->y*n->s[W]->x > 0);
  497.     case Infini3     : return 1;
  498.     case Impossible  : return 1;
  499.     }
  500.     return 2;
  501. }
  502.  
  503.  
  504. int delaunay_tree::conflit(noeud n, DT_item p)
  505. {
  506.   /* teste si un point est en conflit avec une region */
  507.  
  508.     register coord a1,c2,c3,z0;
  509.  
  510.     switch(n->type){
  511.     case Fini       :
  512.         a1 = p->x - n->s[0]->x ;
  513.         c3 = p->y - n->s[0]->y ;
  514.  
  515.         z0 = a1 * (p->x + n->s[0]->x) + c3 * (p->y + n->s[0]->y) ;
  516.  
  517.         c2 = n->a3 * c3 - a1 * n->b3 ;
  518.         c3 = a1 * n->b2 - n->a2 * c3 ;
  519.  
  520.             /********** Arithmetique double precision ************/
  521.             /*      On utilise le fait que les flottants de l'accelerateur
  522.                     flottant sont stockes avec une mantisse de 64 bits !!!
  523.                     Pour faire du calcul entier juste, faite le en flottant !
  524.             */
  525.  
  526.         a1=(double)z0*(double)n->c1 + (double)n->z1*(double)c2 +
  527.                                     (double)n->z2*(double)c3;
  528.         return (n->direct ? ( a1 <= 0.0)  : (a1 >= 0.0) );
  529.     case Infini1    :
  530.         if ( ( (p->x - n->s[U]->x)* (n->s[U]->y - n->s[V]->y)
  531.                 + (p->y - n->s[U]->y)* (n->s[V]->x - n->s[U]->x) ) == 0 )
  532.         {               z0 = n->s[U]->x - p->x;
  533.                         a1 = n->s[U]->y - p->y;
  534.                         c2 = n->s[V]->x - p->x;
  535.                         c3 = n->s[V]->y - p->y;
  536.                         return ( z0*c2 + a1*c3 < 0);
  537.         }else return (n->direct)
  538.             ?   ( ( (p->x - n->s[U]->x)* (n->s[U]->y - n->s[V]->y)
  539.                     + (p->y - n->s[U]->y)* (n->s[V]->x - n->s[U]->x) ) > 0 )
  540.             :   ( ( (p->x - n->s[U]->x)* (n->s[V]->y - n->s[U]->y)
  541.                     + (p->y - n->s[U]->y)* (n->s[U]->x - n->s[V]->x) ) > 0 ) ;
  542.     case Infini2    :
  543.                 return ( (p->x - n->s[U]->x) * (n->s[V]->x + n->s[W]->x)
  544.                     +(p->y - n->s[U]->y) * (n->s[V]->y + n->s[W]->y) > 0 );
  545.     case Infini3    : return 1;
  546.     case Impossible : return 0;
  547.     }
  548.     return 0;
  549. }
  550.  
  551. //int delaunay_tree::conflit( noeud n, DT_item p)
  552. ///* teste si un point est en conflit avec une region */
  553. //{
  554. //    register coord a1,c2,c3,z0;
  555. //
  556. //    switch(n->type){
  557. //    case Fini        :
  558. //        a1 = p->x - n->s[0]->x ;
  559. //        c3 = p->y - n->s[0]->y ;
  560. //
  561. //        z0 = a1 * (p->x + n->s[0]->x) + c3 * (p->y + n->s[0]->y) ;
  562. //     
  563. //        c2 = n->a3 * c3 - a1 * n->b3 ;
  564. //        c3 = a1 * n->b2 - n->a2 * c3 ;
  565. // 
  566. //    /********** Arithmetique double precision ************/
  567. //    /* On utilise le fait que les flottants de l'accelerateur
  568. //       flottant sont stockes avec une mantisse de 64 bits !!!
  569. //           Pour faire du calcul entier juste, faite le en flottant !
  570. //    */
  571. // 
  572. //    a1=(double)z0*(double)n->c1 + (double)n->z1*(double)c2 + (double)n->z2*(double)c3;
  573. //        return (n->direct ? ( a1 <= 0.0)  : (a1 >= 0.0) );
  574. //
  575. //    case Infini1    : 
  576. //        return (n->direct)
  577. //            ?    ( ( (p->x - n->s[U]->x)* (n->s[U]->y - n->s[V]->y)
  578. //                    + (p->y - n->s[U]->y)* (n->s[V]->x - n->s[U]->x) ) >= 0 )
  579. //            :    ( ( (p->x - n->s[U]->x)* (n->s[V]->y - n->s[U]->y)
  580. //                    + (p->y - n->s[U]->y)* (n->s[U]->x - n->s[V]->x) ) >= 0 ) ;
  581. //    case Infini2    :
  582. //                return ( (p->x - n->s[U]->x) * (n->s[V]->x + n->s[W]->x)
  583. //                    +(p->y - n->s[U]->y) * (n->s[V]->y + n->s[W]->y) >= 0 );
  584. //    case Infini3    : return 1;
  585. //    case Impossible : return 0;
  586. //    }
  587. //    return 0;
  588. //}
  589.  
  590.  
  591.  
  592. void delaunay_tree::clear()
  593. {
  594.   if (arbre != nouveau)
  595.   { flag++;
  596.     arbre->s[U]->visite = flag;
  597.     arbre->s[V]->visite = flag;
  598.     arbre->s[W]->visite = flag;
  599.     clear_items(nouveau);
  600.    }
  601.   free_blocks();
  602.   init();
  603. }
  604.  
  605. delaunay_tree::~delaunay_tree() { clear(); }
  606.  
  607. delaunay_tree::delaunay_tree()  { init(); }
  608.  
  609. void delaunay_tree::init()
  610. /* initialisation du Delaunay tree avec les trois points a l'infini */
  611. {
  612.         counter = 0;
  613.         flag=1; 
  614.         Star = NULL; 
  615.         Reinsert = NULL; 
  616.         Poubelle_noeud = Libre_noeud = NULL;
  617.  
  618.         used_blocks=NULL;
  619.         Poubelle_point=NULL;
  620.         Poubelle_listenoeud=NULL;
  621.         Poubelle_noeud=NULL;
  622.         Libre_noeud=NULL;
  623.         Poubelle_star=NULL;
  624.         Poubelle_reinserer=NULL;
  625.         fl = 0;
  626.         item_nombre = 0;
  627.         item_next = 0;
  628.         listenoeud_nombre = 0;
  629.         listenoeud_next = 0;
  630.         noeud_nombre = 0;
  631.         noeud_next = 0;
  632.         star_nombre = 0;
  633.         star_next = 0;
  634.         reinserer_nombre = 0;
  635.         reinserer_next = 0;
  636.         first_star=NULL;
  637.  
  638.  
  639.     arbre = alloc_noeud();
  640.     nouveau = arbre;
  641.  
  642.     arbre->type = Infini3;
  643.     arbre->degenere = 0;
  644.     arbre->visite = 0;
  645.     arbre->stacked = 0;
  646.     arbre->detruire = 0;
  647.     arbre->meurtrier = NULL;
  648.     arbre->fils[U] = arbre->fils[V] = arbre->fils[W] = NULL ;
  649.     arbre->bfils = NULL;
  650.     arbre->pere = NULL;
  651.     arbre->bpere = NULL;
  652.     arbre->voisin[U] =
  653.     arbre->voisin[V] =
  654.     arbre->voisin[W] =
  655.     arbre->creation[U] =
  656.     arbre->creation[V] =
  657.     arbre->creation[W] = alloc_noeud();
  658.     arbre->direct = 1;
  659.  
  660.     arbre->voisin[W]->visite = 0;
  661.     arbre->voisin[W]->stacked = 0;
  662.     arbre->voisin[W]->type = Impossible ;
  663.     arbre->voisin[W]->degenere = 0;
  664.     arbre->voisin[W]->meurtrier = NULL ;
  665.     arbre->voisin[W]->fils[U] = arbre->voisin[W]->fils[V] =
  666.                             arbre->voisin[W]->fils[W] = NULL ;
  667.     arbre->voisin[W]->bfils = NULL ;
  668.     arbre->voisin[W]->pere = NULL ;
  669.     arbre->voisin[W]->bpere = NULL ;
  670.     arbre->voisin[W]->voisin[U] =
  671.     arbre->voisin[W]->voisin[V] =
  672.     arbre->voisin[W]->voisin[W] =
  673.     arbre->voisin[W]->creation[U] =
  674.     arbre->voisin[W]->creation[V] =
  675.     arbre->voisin[W]->creation[W] = arbre ;
  676.     arbre->voisin[W]->direct = 0;
  677.  
  678.     arbre->s[U] = arbre->voisin[W]->s[U] = alloc_point();
  679.     arbre->s[V] = arbre->voisin[W]->s[V] = alloc_point();
  680.     arbre->s[W] = arbre->voisin[W]->s[W] = alloc_point();
  681.     arbre->s[U]->x =  2;
  682.     arbre->s[U]->y =  0;
  683.     arbre->s[V]->x = -1;
  684.     arbre->s[V]->y =  2;
  685.     arbre->s[W]->x = -1;
  686.     arbre->s[W]->y = -2;
  687.     arbre->s[U]->n = 0;
  688.     arbre->s[V]->n = 0;
  689.     arbre->s[W]->n = 0;
  690.  
  691. }
  692.  
  693. noeud delaunay_tree::Insert( noeud n, DT_item p)
  694. /* recherche d'un conflit */
  695. {
  696.     listenoeud l;
  697.     noeud nn;
  698.     int i;
  699.  
  700.     if ( flag != n->visite ) {
  701.      n->visite = flag;
  702.      if ( conflit(n,p)){
  703.         if ( n->meurtrier ) {
  704.             for (i=U; i<=W; i++) if ( n->fils[i] )
  705.                 if ( nn = Insert( n->fils[i], p) ) return nn ;
  706.         } else return n ;
  707.         for ( l = n->bfils; l ; l = l->next )
  708.                 if ( nn = Insert( l->item, p) ) return nn ;
  709.      }
  710.     }
  711.     return NULL ;
  712. }
  713.  
  714. noeud delaunay_tree::creenoeud( noeud pere, int i, DT_item p)
  715. {
  716.                 nouveau = pere->fils[i] = alloc_noeud();
  717.                 beaufils(nouveau, pere->voisin[i] );
  718.                 nouveau->visite = flag;
  719.                 nouveau->stacked = flag;
  720.                 nouveau->detruire = 0;
  721.                 nouveau->meurtrier = NULL;
  722.                 nouveau->fils[U] = nouveau->fils[V] = nouveau->fils[W] = NULL ;
  723.                 nouveau->bfils = NULL;
  724.                 nouveau->s[U]=p;                    /* createur */
  725.                 switch(i){
  726.                 case U: nouveau->direct = pere->direct;
  727.                         nouveau->s[V]= pere->s[V];
  728.                         nouveau->s[W]= pere->s[W]; break;
  729.                 case V: nouveau->direct = ! pere->direct;
  730.                         nouveau->s[V]= pere->s[U];
  731.                         nouveau->s[W]= pere->s[W]; break;
  732.                 case W: nouveau->direct = pere->direct;
  733.                         nouveau->s[V]= pere->s[U];
  734.                         nouveau->s[W]= pere->s[V]; break;
  735.                 }
  736.                 nouveau->pere = pere;
  737.                 nouveau->bpere =
  738.                 nouveau->voisin[U] =
  739.                 nouveau->creation[U] = pere->voisin[i];
  740.                                         /* nouveau est-il fini ou infini */
  741.                 switch ( pere->type ) {
  742.                 case Fini       :
  743.                     cercle(nouveau);
  744.                     nouveau->type = Fini; break;
  745.                 case Infini1    :
  746.                     switch (i) {
  747.                     case U: 
  748.                     case V: demiplan(nouveau); nouveau->type = Infini1;
  749.                             nouveau->degenere = 0; break;
  750.                     case W: cercle(nouveau); nouveau->type = Fini; break;
  751.                     }break;
  752.                 case Infini2    :
  753.                     switch (i) {
  754.                     case U: nouveau->type = Infini2;
  755.                             nouveau->degenere = 0;break; 
  756.                     case V:
  757.                     case W: nouveau->type = Infini1; demiplan(nouveau);
  758.                             nouveau->degenere = 0;break; 
  759.                     }break;
  760.                 case Infini3    : nouveau->type = Infini2;
  761.                                     nouveau->degenere = 0;break; 
  762.                 case Impossible : break;
  763.                 }
  764.     return nouveau;
  765. }
  766.  
  767. int delaunay_tree::creation( noeud detruit, DT_item p)
  768. /* creation eventuelle des nouveaux noeuds et des liens et des voisinages*/
  769. {
  770.     DT_item first= NULL;
  771.     noeud n1,father;
  772.     int arete;        /* indice sur n1 de l'arete axe-createur */ 
  773.     int a,b,c;        /* indice des trois sommets sur father */
  774.     int a1,b1,c1;    /* indice des trois sommets sur le noeud precedent */
  775.     int j;
  776.  
  777.     if ( ! detruit)
  778.                     return 0;
  779.         for (j=U; j<= W - detruit->type; j++) if (confondu(p,detruit->s[j]))
  780.                                     return 0;
  781.         detruit->meurtrier = p;
  782.         a = U; b = V; c = W;
  783.         while ( detruit->voisin[c]->meurtrier
  784.             || conflit ( detruit->voisin[c],p) )
  785.                 {
  786.            detruit->voisin[c]->meurtrier = p;
  787.            a1 = a; b1 = b; c1 = c;
  788.            for (j=U; j<=W; j++)
  789.            if (detruit->voisin[c1]->s[j] == detruit->s[a1]) 
  790.                        a=j;
  791.            else 
  792.                       if (detruit->voisin[c1]->s[j] == detruit->s[b1]) 
  793.                          c=j;
  794.                       else
  795.                          b=j;
  796.            detruit = detruit->voisin[c1] ;
  797.         }
  798.                 p->cree = n1 = creenoeud(detruit,c,p);
  799.                 /* n1 remplace pere comme voisin de son beau-pere*/
  800.                 for (j=U; j<=W; j++)
  801.                     if (    (detruit->voisin[c]->s[j] != n1->s[V])
  802.                         &&    (detruit->voisin[c]->s[j] != n1->s[W])) break;
  803.                 detruit->voisin[c]->voisin[j] = n1 ;
  804.                 switch (c) {
  805.                     case U : first = detruit->s[a=(detruit->direct) ? W : V ];
  806.                                     break;
  807.                     case V : first = detruit->s[a=(detruit->direct) ? U : W ];
  808.                                     break;
  809.                     case W : first = detruit->s[a=(detruit->direct) ? V : U ];
  810.                                     break;
  811.                 }
  812.                 b = c; c = U+V+W -b -a ;
  813.                 father = detruit;
  814.                 arete =  (n1->direct) ? V : W ;
  815.     do {
  816.                 /* ac est l'arete par laquelle n1 est fils de father */
  817.                 /* on veut tourner autour de a en faisant bouger c vers b */
  818.                 /* i e on tourne dans le sens direct */
  819.         while (     father->voisin[c]->meurtrier
  820.                 ||    conflit ( father->voisin[c],p)         ){
  821.             father->voisin[c]->meurtrier = p;
  822.             a1 = a; b1 = b; c1 = c;
  823.             for (j=U; j<=W; j++)
  824.                 if (father->voisin[c1]->s[j] == father->s[a1])            a=j ;
  825.                 else if (father->voisin[c1]->s[j] == father->s[b1])        c=j ;
  826.                 else                                                    b=j;
  827.             father = father->voisin[c1] ;
  828.         }
  829.         if (father->fils[c])
  830.                 n1->creation[arete] = n1->voisin[arete] = father->fils[c] ;
  831.         else {  n1->creation[arete] = n1->voisin[arete] = creenoeud(father,c,p);
  832.                 /* le nouveau remplace pere comme voisin de son beau-pere*/
  833.                 for (j=U; j<=W; j++)
  834.                     if (    (father->voisin[c]->s[j] != father->fils[c]->s[V])
  835.                         &&    (father->voisin[c]->s[j] != father->fils[c]->s[W]))
  836.                                                                         break;
  837.                 father->voisin[c]->voisin[j] = father->fils[c] ;
  838.         }
  839.         arete =  (father->fils[c]->direct) ? V : W ;
  840.         father->fils[c]->creation[V+W-arete] =
  841.                         father->fils[c]->voisin[V+W-arete] = n1 ;
  842.         n1 = father->fils[c] ;
  843.         j = a ;
  844.         a = b ;
  845.         b = c ;
  846.         c = j ;
  847.     } while ( father->s[a] != first );
  848.     return 1;
  849. }
  850.  
  851. DT_item delaunay_tree::localise( noeud detecte, DT_item p)
  852. /* determination du plus proche voisin */
  853. {
  854.     DT_item first= NULL,trouve = NULL;
  855.     coord dist, distance = 10000000.0;
  856.     noeud father;
  857.     int a,b,c;        /* indice des trois sommets sur father */
  858.     int a1,b1,c1;    /* indice des trois sommets sur le noeud precedent */
  859.     int j;
  860.  
  861.     if (!detecte)
  862.                     return NULL;
  863.  
  864.         for (j=U; j<= W - detecte->type; j++) if (confondu(p,detecte->s[j]))
  865.                                     return detecte->s[j];
  866.         detecte->visite = flag++ ;
  867.         a = U; b = V; c = W;
  868.         while (     ( detecte->voisin[c]->visite == flag )
  869.                 ||    conflit ( detecte->voisin[c],p)         ){
  870.             detecte->voisin[c]->visite = flag;
  871.             a1 = a; b1 = b; c1 = c;
  872.             for (j=U; j<=W; j++)
  873.                 if (detecte->voisin[c1]->s[j] == detecte->s[a1])        a=j ;
  874.                 else if (detecte->voisin[c1]->s[j] == detecte->s[b1])    c=j ;
  875.                 else                                                    b=j;
  876.             detecte = detecte->voisin[c1] ;
  877.         }
  878.                 switch (c) {
  879.                     case U : first = detecte->s[a=(detecte->direct) ? W : V ];
  880.                                     break;
  881.                     case V : first = detecte->s[a=(detecte->direct) ? U : W ];
  882.                                     break;
  883.                     case W : first = detecte->s[a=(detecte->direct) ? V : U ];
  884.                                     break;
  885.                 }
  886.                 b = c; c = U+V+W -b -a ;
  887.                 father = detecte;
  888.     do {
  889.                 /* on veut tourner autour de a en faisant bouger c vers b */
  890.                 /* i e on tourne dans le sens direct */
  891.         while (     ( father->voisin[c]->visite == flag )
  892.                 ||    conflit ( father->voisin[c],p)         ){
  893.             father->voisin[c]->visite = flag;
  894.             a1 = a; b1 = b; c1 = c;
  895.             for (j=U; j<=W; j++)
  896.                 if (father->voisin[c1]->s[j] == father->s[a1])            a=j ;
  897.                 else if (father->voisin[c1]->s[j] == father->s[b1])        c=j ;
  898.                 else                                                    b=j;
  899.             father = father->voisin[c1] ;
  900.         }
  901.         if (father->s[a]->n) {
  902.             dist = (p->x - father->s[a]->x) * (p->x - father->s[a]->x) +
  903.                     (p->y - father->s[a]->y) * (p->y - father->s[a]->y)    ;
  904.             if (dist < distance) { distance = dist; trouve = father->s[a]; }
  905.         }
  906.         j = a ;
  907.         a = b ;
  908.         b = c ;
  909.         c = j ;
  910.     } while ( father->s[a] != first );
  911.     return trouve;
  912. }
  913.  
  914. void delaunay_tree::add_x1x2( reinserer r, noeud v, DT_item p)
  915. /* ajoute v a comme triangle detruit de r */
  916. {
  917.     if ( v->direct )
  918.                      if  (p==v->s[V]) r->detruit1 = v;
  919.                      else              r->detruit2 = v;
  920.     else            
  921.                      if  (p==v->s[W]) r->detruit1 = v;
  922.                      else              r->detruit2 = v;
  923. }
  924.  
  925. void delaunay_tree::add_decroche( reinserer r, noeud v)
  926. /* ajoute v a la liste des decroches de r */
  927. {
  928.     listenoeud nov;
  929.  
  930.     nov = alloc_listenoeud();
  931.     nov->next = r->decroches;
  932.     nov->item = v;
  933.     r->decroches = nov;
  934. }
  935.  
  936. void delaunay_tree::destroy( noeud u, DT_item p)
  937. /* recherche recursive pour initialisation de Reinsert */
  938. {
  939.     int i;
  940.     listenoeud l;
  941.     reinserer rr;
  942.  
  943.     if ( ! u ) return;
  944.     if (u->detruire == flag) return;
  945.     for ( i=U ; i<=W; i++ ) if (u->s[i]==p) break;
  946.     if ( i<=W ) {
  947.         u->detruire = flag;
  948.         if (u->meurtrier)
  949.             for ( i=U ; i<=W; i++ ) if (u->fils[i]) destroy(u->fils[i], p);
  950.         for ( l=u->bfils; l; l = l->next )      destroy(l->item   , p);
  951.         if (u->s[U]==p) return;
  952.         rr = locate(Reinsert,u->s[U]);
  953.         add_x1x2(rr,u,p);
  954.     } else {
  955.         rr = locate(Reinsert,u->s[U]);
  956.         add_decroche(rr,u);
  957.         if (u->pere->detruire == flag)    u->pere        = NULL;
  958.         else                            u->bpere    = NULL;
  959.     }
  960. }
  961.  
  962. void delaunay_tree::recherche( DT_item p)
  963. /* init des structures Star et Reinsert */
  964. {
  965.     DT_item first;
  966.     noeud father;
  967.     int arete;        /* indice sur n1 de l'arete axe-createur */ 
  968.     int a,b,c;        /* indice des trois sommets sur father */
  969.     int a1,b1,c1;    /* indice des trois sommets sur le noeud precedent */
  970.     int j;
  971.  
  972.     clear_Reinsert(Reinsert); 
  973.     clear_Star();
  974.         arete = (  p->cree->direct  ) ? V : W ;
  975.         first = p->cree->s[V + W - arete ];
  976.         father = p->cree->pere;
  977.         for (j=U; j<=W; j++)
  978.             if (father->s[j] == p->cree->s[V + W - arete ])        a=j ;
  979.             else if (father->s[j] == p->cree->s[arete])            c=j ;
  980.             else                                b=j;
  981.     do {
  982.                 /* ac est l'arete par laquelle on est fils de father */
  983.                 /* on veut tourner autour de a en faisant bouger c vers b */
  984.                 /* et que on tourne autour de a dans le sens inverse */
  985.         while ( ( father->voisin[c]->meurtrier == p )
  986.                 || (father->voisin[c]->visite == flag) ) {
  987.             a1 = a; b1 = b; c1 = c;
  988.             for (j=U; j<=W; j++)
  989.                 if (father->voisin[c1]->s[j] == father->s[a1])            a=j ;
  990.                 else if (father->voisin[c1]->s[j] == father->s[b1])        c=j ;
  991.                 else                                                    b=j;
  992.             father = father->voisin[c1] ;
  993.             father->visite = flag ;
  994.             father->meurtrier = NULL;
  995.         }
  996.         father->visite = flag ;
  997.         father->meurtrier = NULL;
  998.         a1 = a; b1 = b; c1 = c;
  999.         for (j=U; j<=W; j++)
  1000.                 if (father->voisin[c1]->s[j] == father->s[a1])            a=j ;
  1001.                 else if (father->voisin[c1]->s[j] == father->s[b1])        c=j ;
  1002.                 else                                                    b=j;
  1003.         father->voisin[c1]->special[b] = father;
  1004.         if ( father->voisin[c1]->voisin[b]->s[U] == p )
  1005.                                 father->voisin[c1]->voisin[b] = father;
  1006.         plusbeaufils(father->fils[c1], father->voisin[c1] );
  1007.         destroy(father->fils[c1],p);
  1008.         poubelle_noeud(father->fils[c1]);
  1009.         father->fils[c1] = NULL;
  1010.         init_Star(father,a1,b1,c1);
  1011.         a = b1 ;
  1012.         b = c1 ;
  1013.         c = a1 ;
  1014.     } while ( father->s[a] != first );
  1015.     init_Star( (noeud) NULL, 0, 0, 0);
  1016. }
  1017.  
  1018. void delaunay_tree::reinsertion( reinserer n, DT_item p)
  1019. {
  1020.     listenoeud l;
  1021.     noeud u,v,w,ww,w1,w2;
  1022.     DT_item x1,x2;
  1023.     int i,j,k,a,b,c,a1,b1,c1;
  1024.     star s1,s2;
  1025.  
  1026.     /* traitement des noeuds decroches */
  1027.     for ( l=n->decroches; l; l = l->next )
  1028.         if ( l->item->pere ) {
  1029.         /* il faut trouver un nouveau beau-pere et tous les voisinages 
  1030.                    par cette arete */
  1031.             u = l->item->pere;
  1032.             for( i=U; i<=W; i++) if (u->fils[i] == l->item) break;
  1033.             l->item->bpere = u->special[i];
  1034.             l->item->creation[U] = u->special[i];
  1035.             l->item->special[U] = u->special[i];
  1036.             if ( l->item->voisin[U]->detruire == flag )
  1037.                                 l->item->voisin[U] = u->special[i];
  1038.             beaufils(l->item, u->special[i]);
  1039.             for( j=U; j<=W; j++)
  1040.                 if (    (l->item->bpere->s[j] != l->item->s[V] )
  1041.                     &&    (l->item->bpere->s[j] != l->item->s[W] ) ) break;
  1042.             l->item->bpere->voisin[j] = l->item;
  1043.         } else {
  1044.         /* il faut trouver un nouveau pere */
  1045.             u = l->item->bpere;
  1046.             for( i=U; i<=W; i++) if (u->s[i] == l->item->s[V]) break;
  1047.             for( j=U; j<=W; j++) if (u->s[j] == l->item->s[W]) break;
  1048.             i = U+V+W -i -j;
  1049.             u= l->item->pere = u->special[i];
  1050.             for( i=U; i<=W; i++) if (u->s[i] == l->item->s[V]) break;
  1051.             for( j=U; j<=W; j++) if (u->s[j] == l->item->s[W]) break;
  1052.             i = U+V+W -i -j;
  1053.             u->fils[i] = l->item;
  1054.         }
  1055.     /* traitement des noeuds detruits */
  1056.     if ( ! n->detruit1 ) return;
  1057.     x1 =  ( n->detruit1->s[V] == p ) ? n->detruit1->s[W] : n->detruit1->s[V];
  1058.     x2 =  ( n->detruit2->s[V] == p ) ? n->detruit2->s[W] : n->detruit2->s[V];
  1059.     s1 = loc_Star(x1);
  1060.     s2 = loc_Star(x2)->previous ;
  1061.     u = s1->triangle;
  1062.     /* Cas 1 (triangle x x1 x2 a creer) u n'est pas en conflit avec x */
  1063.     if ( ! conflit(u, n->x) ) {
  1064.             /* u est le beau pere de x x1 x2 a cree */
  1065.             for( i=U; i<=W; i++) if (u->s[i] == x1) break;
  1066.             for( j=U; j<=W; j++) if (u->s[j] == x2) break;
  1067.             k = U+V+W -i -j;
  1068.             u = u->voisin[k] ;
  1069.             /* u est le pere */
  1070.             for( i=U; i<=W; i++) if (u->s[i] == x1) break;
  1071.             for( j=U; j<=W; j++) if (u->s[j] == x2) break;
  1072.             i = U+V+W -i -j;
  1073.             n->x->cree =ww = creenoeud(u ,i,n->x);
  1074.         /* ww est le fils de son pere */
  1075.             u->voisin[i] = ww->bpere;        /* voisin a la mort */
  1076.         /* voisinage beau-fils beau-pere */
  1077.             ww->bpere->voisin[k] = ww;
  1078.         /* voisinage avec les "freres" du nouveau */
  1079.             i = ( ww->s[V] == x1 ) ? V : W ;
  1080.             u = n->detruit2->creation[(n->detruit2->s[V]==p)?V:W] ;
  1081.             ww->voisin[i] = ww->creation[i] = u;
  1082.             for( j=U; j<=W; j++) if ( u->creation[j] == n->detruit2 ) break;
  1083.             u->creation[j] = u->special[j] = ww ;
  1084.             if ( u->voisin[j]->detruire == flag ) u->voisin[j] = ww ;
  1085.             if ( u->voisin[j]->visite == flag ) u->voisin[j] = ww ;
  1086.             u = n->detruit1->creation[(n->detruit1->s[V]==p)?V:W] ;
  1087.             ww->voisin[V+W-i] = ww->creation[V+W-i] = u;
  1088.             for( j=U; j<=W; j++) if ( u->creation[j] == n->detruit1 ) break;
  1089.             u->creation[j] = u->special[j] = ww ;
  1090.             if ( u->voisin[j]->detruire == flag ) u->voisin[j] = ww ;
  1091.             if ( u->voisin[j]->visite == flag ) u->voisin[j] = ww ;
  1092.         /* mise a jour de Star */
  1093.             poubelle_noeud(n->detruit1);
  1094.             poubelle_noeud(n->detruit2);
  1095.             split_Star(s1,s2); s2=s1->next;
  1096.             s1->triangle = ww;
  1097.             s2->triangle = ww;
  1098.             for( i=U; i<=W; i++) if (ww->s[i] == x1)    s1->p = s2->arete = i;
  1099.                         else     if (ww->s[i] == n->x ) s1->n = s2->p = i;
  1100.                         else                           s1->arete = s2->n = i;
  1101.             ww->Star[s1->arete] = s1;
  1102.             ww->Star[s2->arete] = s2;
  1103.     return;
  1104.     }
  1105.     /* Cas 2  */
  1106.     v = u; w1=NULL;
  1107.     w = n->detruit1->creation[(n->detruit1->s[V]==p)?V:W] ;
  1108.     a = s1->p ;
  1109.     c = s1->n ;
  1110.     b = s1->arete ;
  1111.     do {
  1112.     /* on tourne dans le sens inverse */
  1113.         while ( conflit( v->voisin[c], n->x ) ) {
  1114.             a1 = a; b1 = b; c1 = c;
  1115.             for (j=U; j<=W; j++) if ( v->voisin[c1]->s[j] == v->s[a1]) a = j;
  1116.                         else     if ( v->voisin[c1]->s[j] == v->s[b1]) c = j;
  1117.                         else                                           b = j;
  1118.             v = v->voisin[c1];
  1119.             v->meurtrier = n->x;
  1120.         }
  1121.             v->meurtrier = n->x;
  1122.             ww = creenoeud(v,c,n->x);
  1123.         /* voisinage beau-fils beau-pere */
  1124.             for (j=U; j<=W; j++)
  1125.                     if (    ( ww->bpere->s[j] != ww->s[V] )
  1126.                         &&    ( ww->bpere->s[j] != ww->s[W] ) ) break;
  1127.             if ( ww->bpere->visite == flag ){ /* le beau-pere est dans A */
  1128.                 ww->bpere->voisin[j] = ww;
  1129.             } else {
  1130.                 ww->bpere->special[j] = ww;
  1131.                 if ( ww->bpere->voisin[j]->detruire == flag )
  1132.                                                 ww->bpere->voisin[j] = ww ;
  1133.                 if ( ww->bpere->voisin[j]->visite == flag )
  1134.                                                 ww->bpere->voisin[j] = ww ;
  1135.             }
  1136.         /* voisinage nouvelle arete */
  1137.             for (j=U;j<=W;j++) if((ww->s[j]!=n->x)&&(ww->s[j]!=v->s[a]))break;
  1138.             ww->voisin[j] = ww->creation[j] = w ;
  1139.             if ( !w1 ) {
  1140.               for (j=U;j<=W;j++) if((w->s[j]!=n->x)&&(w->s[j]!=v->s[a]))break;
  1141.               w1 = w->special[j] = w->creation[j] = ww ;
  1142.               if (w->voisin[j]->visite == flag ) w->voisin[j] = ww ;
  1143.               if (w->voisin[j]->detruire == flag ) w->voisin[j] = ww ;
  1144.             } else {
  1145.               for (j=U;j<=W;j++) if((w->s[j]!=n->x)&&(w->s[j]!=v->s[a]))break;
  1146.               w->voisin[j] = w->creation[j] = ww ;
  1147.             }
  1148.             w = ww;
  1149.         /* mise a jour de Star */
  1150.             if ( ww->bpere->visite != flag ){/* le beau-pere n'est pas dans A */
  1151.                 v->Star[c]->triangle = ww;
  1152.                 ww->Star[U] = v->Star[c];
  1153.                 ww->Star[U]->arete = U;
  1154.                 ww->Star[U]->n = (v->s[a] == w->s[V] ) ? V : W ;
  1155.                 ww->Star[U]->p = V+W-ww->Star[U]->n ;
  1156.             }
  1157.         j = a ; a = b ; b = c ; c = j ;
  1158.     } while ( v->s[a] != x2) ; 
  1159.         /* voisinage autour de x x2 */
  1160.             n->x->cree = w2 = w;
  1161.             ww = n->detruit2->creation[(n->detruit2->s[V]==p)?V:W] ;
  1162.             for (j=U; j<=W; j++) if((ww->s[j]!=n->x)&&(ww->s[j]!=v->s[a]))break;
  1163.             if (ww->voisin[j]->visite == flag ) ww->voisin[j] = w ;
  1164.             if (ww->voisin[j]->detruire == flag ) ww->voisin[j] = w ;
  1165.             ww->creation[j] = ww->special[j] = w ;
  1166.             for (j=U; j<=W; j++) if((w->s[j]!=n->x)&&(w->s[j]!=v->s[a]))break;
  1167.             w->voisin[j] = w->creation[j] = ww ;
  1168.     do {
  1169.     /* on tourne toujours dans le sens inverse */
  1170.         while ( conflit( v->voisin[c], n->x ) ) {
  1171.             a1 = a; b1 = b; c1 = c;
  1172.             for (j=U; j<=W; j++) if ( v->voisin[c1]->s[j] == v->s[a1]) a = j;
  1173.                         else     if ( v->voisin[c1]->s[j] == v->s[b1]) c = j;
  1174.                         else                                           b = j;
  1175.             v = v->voisin[c1];
  1176.             v->meurtrier = n->x;
  1177.         }
  1178.             v->meurtrier = n->x;
  1179.         j = a ; a = b ; b = c ; c = j ;
  1180.     } while ( v->s[a] != x1) ; 
  1181.         /* mise a jour de Star */
  1182.             split_Star(s1,s2); s2=s1->next;
  1183.             s1->triangle = w1;
  1184.             s2->triangle = w2;
  1185.             poubelle_noeud(n->detruit1);
  1186.             poubelle_noeud(n->detruit2);
  1187.             for( i=U; i<=W; i++) if (w1->s[i] == x1)    s1->p = i;
  1188.                         else     if (w1->s[i] == n->x ) s1->n = i;
  1189.                         else                            s1->arete = i;
  1190.             w1->Star[s1->arete] = s1;
  1191.             for( i=U; i<=W; i++) if (w2->s[i] == x2)    s2->n = i;
  1192.                         else     if (w2->s[i] == n->x ) s2->p = i;
  1193.                         else                            s2->arete = i;
  1194.             w2->Star[s2->arete] = s2;
  1195. }
  1196.  
  1197. void delaunay_tree::parcours( reinserer n, DT_item p)
  1198. {
  1199.     if ( ! n ) return;
  1200.     parcours(n->finf,p);
  1201.     reinsertion(n,p);
  1202.     parcours(n->fsup,p);
  1203. }
  1204.  
  1205.  
  1206. void delaunay_tree::trace_items(noeud n, delaunay_f2* trace_item)
  1207. /* exploration recursive pour le trace des points */
  1208. {
  1209.   noeud* stack = new noeud[counter+10];
  1210.   int top = 0;
  1211.   stack[0] = n;
  1212.   n->stacked = flag;
  1213.  
  1214.   while (top >= 0)
  1215.   { n = stack[top--];   //pop
  1216.     n->visite = flag;
  1217.     for (int i=U; i<=W; i++)
  1218.       if ( n->s[i] && n->s[i]->visite != flag )
  1219.             { n->s[i]->visite = flag ;
  1220.               trace_item(n->s[i]->x, n->s[i]->y);
  1221.              }
  1222.  
  1223.     for (i=W; i>=U; i--)
  1224.         { noeud v = n->voisin[i];
  1225.       if ( v->stacked != flag )
  1226.            { stack[++top] = v;          //push
  1227.              v->stacked = flag;
  1228.             }
  1229.          }
  1230.  
  1231.    } // while stack not empty
  1232.  
  1233.  delete stack;
  1234.  
  1235. }
  1236.  
  1237. void delaunay_tree::list_items(noeud n, list<DT_item>& L)
  1238. /* exploration recursive pour le trace des points */
  1239. {
  1240.   noeud* stack = new noeud[counter+10];
  1241.   int top = 0;
  1242.   stack[0] = n;
  1243.   n->stacked = flag;
  1244.  
  1245.   while (top >= 0)
  1246.   { n = stack[top--];   //pop
  1247.     n->visite = flag;
  1248.     for (int i=U; i<=W; i++)
  1249.       if ( n->s[i] && n->s[i]->visite != flag )
  1250.             { n->s[i]->visite = flag ;
  1251.               L.append(n->s[i]);
  1252.              }
  1253.  
  1254.     for (i=W; i>=U; i--)
  1255.         { noeud v = n->voisin[i];
  1256.       if ( v->stacked != flag )
  1257.            { stack[++top] = v;          //push
  1258.              v->stacked = flag;
  1259.             }
  1260.          }
  1261.  
  1262.    } // while stack not empty
  1263.  
  1264.  delete stack;
  1265.  
  1266. }
  1267.  
  1268. void delaunay_tree::clear_items(noeud n)
  1269. {
  1270.   noeud* stack = new noeud[counter+10];
  1271.   int top = 0;
  1272.   stack[0] = n;
  1273.   n->stacked = flag;
  1274.  
  1275.   while (top >= 0)
  1276.   { n = stack[top--];   //pop
  1277.     n->visite = flag;
  1278.     for (int i=U; i<=W; i++)
  1279.       if ( n->s[i] && n->s[i]->visite != flag )
  1280.             { n->s[i]->visite = flag ;
  1281.               clear_inf(n->s[i]->i);
  1282.              }
  1283.  
  1284.     for (i=W; i>=U; i--)
  1285.         { noeud v = n->voisin[i];
  1286.       if ( v->stacked != flag )
  1287.            { stack[++top] = v;          //push
  1288.              v->stacked = flag;
  1289.             }
  1290.          }
  1291.  
  1292.    } // while stack not empty
  1293.  
  1294.  delete stack;
  1295. }
  1296.  
  1297.  
  1298. void delaunay_tree::del_noeud( noeud n, delaunay_f4* trace_segment,int both)
  1299. /* exploration recursive de delaunay_tree */
  1300. {
  1301.   noeud* stack = new noeud[counter+10];
  1302.   int top = 0;
  1303.   stack[0] = n;
  1304.   n->stacked = flag;
  1305.  
  1306.   int i;
  1307.  
  1308.   while (top >= 0)
  1309.   { n = stack[top--];   //pop
  1310.  
  1311.     n->visite = flag;
  1312.  
  1313.     for (i=U; i<=W; i++)
  1314.           if ( both || n->voisin[i]->visite != flag ) 
  1315.            { int j = (i==U) ? V : U ; 
  1316.              int k = U + V + W -i -j ;
  1317.          if (( n->type == Fini ) || ((i == W) && (n->type==Infini1) ) )
  1318.                 trace_segment(n->s[j]->x,n->s[j]->y, n->s[k]->x, n->s[k]->y);
  1319.         }
  1320.  
  1321.     for (i=W; i>=U; i--)
  1322.         { noeud v = n->voisin[i];
  1323.       if ( v->stacked != flag )
  1324.            { stack[++top] = v;          //push
  1325.              v->stacked = flag;
  1326.             }
  1327.          }
  1328.  
  1329.    } // while stack not empty
  1330.  
  1331.  delete stack;
  1332.  
  1333.  }
  1334.  
  1335.  
  1336. void delaunay_tree::vor( noeud n, delaunay_f6* trace_segment, delaunay_F6* pt_infi, int both)
  1337. {
  1338.  
  1339.   int i,j;
  1340.   coord x,y,mx,my;
  1341.   int order[3];
  1342.  
  1343.   int top = 0;
  1344.   noeud* stack = new noeud[counter+10];
  1345.   stack[0] = n;
  1346.   n->stacked = flag;
  1347.  
  1348.   while (top >= 0)
  1349.   { 
  1350.     n = stack[top--];   //pop
  1351.  
  1352.     n->visite = flag;
  1353.  
  1354.         if (n->direct)
  1355.          { order[0] = 0;
  1356.            order[1] = 1;
  1357.            order[2] = 2;
  1358.           }
  1359.         else 
  1360.          { order[0] = 2;
  1361.            order[1] = 1;
  1362.            order[2] = 0;
  1363.           }
  1364.  
  1365.     for (j=0; j<3; j++)
  1366.         { 
  1367.           int i = order[j];
  1368.  
  1369.                if (both ||  n->voisin[i]->visite != flag )
  1370.         {   int k;
  1371.                     if (n->direct)
  1372.                        k = (i==2) ? 0 : i+1; 
  1373.                     else 
  1374.                        k = (i==0) ? 2 : i-1; 
  1375.  
  1376.             coord sx0 = n->s[k]->x;
  1377.                     coord sy0 = n->s[k]->y;
  1378.  
  1379.             switch (n->type){
  1380.             case Fini :
  1381.             if (! n->degenere ) {
  1382.                 switch ( n->voisin[i]->type ){
  1383.                 case Fini : if ( ! n->voisin[i]->degenere ){
  1384. trace_segment(n->cx,n->cy, n->voisin[i]->cx, n->voisin[i]->cy,sx0,sy0); 
  1385.                                             break;
  1386.                 }
  1387.                 case Infini1  : 
  1388.                   pt_infi(n->cx,n->cy, n->voisin[i]->cx, n->voisin[i]->cy,&x,&y);
  1389.                                   trace_segment(n->cx, n->cy,x,y,sx0,sy0);
  1390.                   break ;
  1391.                 }
  1392.                 break;
  1393.             }
  1394.             case Infini1  :
  1395.                 if ( n->degenere &&  n->voisin[i]->degenere ) break;
  1396.                 switch ( n->voisin[i]->type ){
  1397.                 case Fini :
  1398.                 if ( ! n->voisin[i]->degenere ) 
  1399.                                 { pt_infi(n->voisin[i]->cx, 
  1400.                                           n->voisin[i]->cy, 
  1401.                                           n->cx,n->cy,&x,&y);
  1402.                                   trace_segment(x,y,
  1403.                                                 n->voisin[i]->cx, 
  1404.                                                 n->voisin[i]->cy,
  1405.                                                 sx0,sy0);
  1406.                   break ;
  1407.                  }
  1408.                 case Infini1  :
  1409.                   if (i == W)
  1410.                                   { mx = ( n->s[U]->x + n->s[V]->x )/2 ;
  1411.                     my = ( n->s[U]->y + n->s[V]->y )/2 ;
  1412.                     pt_infi(mx,my, n->cx,n->cy,&x,&y);
  1413.                     pt_infi(mx,my, n->voisin[i]->cx, 
  1414.                                                    n->voisin[i]->cy,&mx, &my);
  1415.                     trace_segment(mx,my, x,y,sx0,sy0);
  1416.                    }
  1417.                 }
  1418.                 break;
  1419.  
  1420.             } //switch
  1421.  
  1422.                   } //if
  1423.  
  1424.            }  // for
  1425.  
  1426.     for (i=W; i>=U; i--)
  1427.         { noeud v = n->voisin[i];
  1428.       if ( v->stacked != flag )
  1429.            { stack[++top] = v;          //push
  1430.              v->stacked = flag;
  1431.             }
  1432.          }
  1433.  
  1434.    } // while stack not empty
  1435.  
  1436.  delete stack;
  1437.  
  1438. }
  1439.  
  1440. DT_item delaunay_tree::neighbor(double x,double y)
  1441. /* recherche du plus proche voisin de (x,y) */
  1442. {
  1443.     DT_POINT requete;
  1444.     DT_item reponse;
  1445.  
  1446.     flag++;
  1447.     if (arbre == nouveau) return NULL;
  1448.     requete.x = x; requete.y = y;
  1449.     reponse = localise(Insert(arbre,&requete), &requete);
  1450.     return  reponse;
  1451. }
  1452.  
  1453. void delaunay_tree::del(double x, double y)
  1454. { DT_item p = neighbor(x,y);
  1455.   if (p==nil || p->x != x || p->y != y) return;
  1456.   del_item(p);
  1457.  }
  1458.  
  1459. void delaunay_tree::del_item( DT_item ppp)
  1460. {
  1461.     DT_item p;
  1462.  
  1463.     if (!ppp) return;
  1464.  
  1465.     flag++;
  1466.     p = (DT_item) ppp ;
  1467.     nouveau = p->cree->pere;    /* pour si on supprime le dernier */
  1468.     recherche(p);
  1469.     last = nouveau ;
  1470.     parcours(Reinsert,p);
  1471.         clear_inf(p->i);
  1472.         counter--;
  1473.     poubelle_point(p);
  1474. }
  1475.  
  1476. DT_item delaunay_tree::insert( coord x, coord y, void* inf)
  1477. /* procedure d'insertion d'un point dans l'arbre */
  1478. {
  1479.     DT_item p;
  1480.  
  1481.     flag++;
  1482.     p = alloc_point();
  1483.     p->x = x;
  1484.     p->y = y;
  1485.         p->i = inf;
  1486.     p->n = flag;
  1487.     if (! creation( Insert(arbre,p), p )) 
  1488.         { poubelle_point(p); 
  1489.           return neighbor(x,y);
  1490.          }
  1491.         copy_inf(inf);
  1492.         counter++;
  1493.     last = nouveau ;
  1494.     return p;
  1495. }
  1496.  
  1497. void delaunay_tree::all_items(list<DT_item>& L)
  1498. /* procedure tracant tous les points */
  1499. {
  1500.     if (arbre == nouveau) return;
  1501.     flag++;
  1502.     arbre->s[U]->visite = flag;
  1503.     arbre->s[V]->visite = flag;
  1504.     arbre->s[W]->visite = flag;
  1505.     list_items(nouveau, L);
  1506. }
  1507.  
  1508. void delaunay_tree::trace_triang_edges( delaunay_f4 *seg,int both)
  1509. /* procedure de trace de delaunay */
  1510. {
  1511.     if (arbre == nouveau) return;
  1512.     flag++;
  1513.     del_noeud(nouveau,seg, both);
  1514. }
  1515.  
  1516. void delaunay_tree::trace_voronoi_sites(delaunay_f2 *trace_item)
  1517. {
  1518.     if (arbre == nouveau) return;
  1519.     flag++;
  1520.     trace_items(nouveau,trace_item);
  1521. }
  1522.  
  1523. void delaunay_tree::trace_voronoi_edges(delaunay_f6 *trace_seg, delaunay_F6 *pt_infi, int both)
  1524. {
  1525.     if (arbre == nouveau) return;
  1526.     flag++;
  1527.     vor(nouveau,trace_seg,pt_infi, both);
  1528. }
  1529.  
  1530. void delaunay_tree::convex_hull(list<DT_item>& CH)
  1531. { 
  1532.  CH.clear();
  1533.  
  1534.  if (arbre == nouveau) return;
  1535.  
  1536.  for (int j=2;j>=0;j--)
  1537.  { noeud v = arbre->voisin[0];
  1538.  
  1539.    int i = (j == 0) ? 2 : j-1;
  1540.  
  1541.    DT_item start = v->s[i];
  1542.  
  1543.    do
  1544.    { int k;
  1545.      if (v->type == Infini1) CH.append(v->s[i]);
  1546.      if (v->direct)
  1547.         k = (i == 0) ? 2 : i-1;
  1548.      else
  1549.         k = (i == 2) ? 0 : i+1;
  1550.      noeud w = v->voisin[k];
  1551.      for(i=0;w->voisin[i] != v;i++); 
  1552.      v = w;
  1553.     } while (v->s[i] != start);
  1554.  
  1555.   }
  1556. }
  1557.  
  1558.